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This article reviews and presents various solved and open problems in the development, 
analysis, and control of epidemic models. Proper modeling and analysis of spreading processes 
has been a longstanding area of research among many different fields including mathematical 
biology, physics, computer science, engineering, economics, and the social sciences. One of 
the earliest epidemic models conceived was by Daniel Bernoulli in 1760, which was motivated 
by studying the spreading of smallpox [Tj. In addition to Bernoulli, there were many different 
researchers also working on mathematical epidemic models around this time [2j. These initial 
models were quite simplistic and the further development and study of such models dates back 
to the 1900s [3|-[j6j, where still simple models were studied to provide insight as to how various 
diseases can spread through a population. In recent years, we have seen a resurgence of interest 
in these problems as the concept of ‘networks’ becomes increasingly prevalent in modeling 
many different aspects of our world today. A more comprehensive review of the history of 
mathematical epidemiology can be found in (7), (8j. 

Despite the study of epidemic models having spanned such a long period of time, it is 
only recently that control engineers have entered the scene. Consequently, there is already a vast 
body of work dedicated to the development and analysis of epidemic models, but far less that 
provide proper insight and machinery on how to effectively control these processes. The focus 
of this article is to provide an introductory tutorial on the latter. We are interested in presenting 
a relatively concise report for new engineers looking to enter the field of spreading processes 
on complex networks. This article presents some of the more well-known and recent results in 
the literature while also identifying numerous open problems that can benefit from the collective 
knowledge of optimization and control theorists. 
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Although this article focuses on the context of epidemics, the same models and tools we 
present are directly applicable to a myriad of different spreading processes on complex networks. 
Examples include the adoption of an idea or rumor through a social network like Twitter, the 
consumption of a new product in a marketplace, the risk of receiving a computer virus through 
the World Wide Web, and of course the spreading of a disease through a population [j9j-[ 11]. For 
this reason, we will often use the terms individuals, people, nodes, and agents interchangeably. 


We begin the article by introducing and analyzing some classical epidemic models and 
their extensions to network models. We then discuss various methods of controlling these 
epidemic models and several extensions. After describing the main shortcomings in the current 
literature and highlighting some recent preliminary works that are aimed at improving the current 
state of the art, we close by providing some intuition into the current research challenges that 
need to be addressed in order to fully harness the power of these works and make a real societal 
impact. 


Modeling and Analysis of Epidemics 


Before jumping into the class of models we study in this article, we must start by 
emphasizing that there are an uncountable number of ways to model spreading processes. 
The underlying common factor that ties almost all epidemic models together is the existence 
of ‘compartments’ in which individuals in a population are divided. The two most common 
compartments that exist in essentially every single epidemic model are called ‘Susceptible’ and 
‘Infected’ (6j, |7J, [1_2|. In models that contain only these two compartments, a given population 
is initially divided into these two compartments. The ‘Susceptible’ compartment (S') represents 
individuals who are healthy, but susceptible to becoming infected. The ’Infected’ compartment (/) 
captures individuals who are infected, but able to recover. From here there are an insurmountable 
number of ways that the interactions within the population can be modeled. 


Throughout this article, we focus on models where individuals can move from one 
compartment to another randomly with some defined rates. For instance, in this two compartment 
model healthy individuals can randomly transition from S to I with some infection rate that is 
a result of interactions with infected individuals. Similarly, infected individuals can randomly 
transition from I to S with some recovery rate that is a result of recovering from the infection. 
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More details on how these rates are defined are provided later. Figure [T] shows the simple 
interaction described above. 


In addition to models with only two compartments, there are many other epidemic models 
aimed at capturing various important features of realistic diseases and spreading processes. This is 
often done by adding more compartments, such as a ‘Removed’ ( R ) compartment representing 
individuals who are no longer susceptible to the infection. This might refer to a deceased, 
vaccinated, or immune individual. Other compartments have also been proposed in the literature 
to study the effect of, for example, an incubation period, partial immunity, or quarantine in the 
spreading dynamics [ 13J—[[T9|. 

For brevity, we will now focus our attention on two of the oldest epidemic models 
known as the Susceptible-Infected-Removed (SIR) and Susceptible-Infected-Susceptible (SIS) 
models 0 . However, we note that the following exposition can be applied to more general 
compartmental models with some appropriate changes. 


Let N be the total number of individuals in a population. We denote by X ^(f) e {S, /, R} 
the state of node i e (1,..., iV} at time t. We collect the state of the entire population in a state 
vector X(t) = (Xi(t),... ,X N (t)) 7 . The evolution of the states is then described by a Markov 
process as follows. An individual i infected at time t recovers to the removed state at a fixed 
rate Si > 0. In other words, if node i is infected at time t, the probability that this node loses 
its infection in the time slot (t , t + At] for small At is given by 

Pr(3Q(t + At) = R\ Xi(t) = I) = SAt + o(At), (SIR) 

Pr (Xi(t + At) = S | Xi(t) = I) = 5 t At + o(At). (SIS) 


The above represents an endogenous transition which occurs internally within each node, 


independent of the states of other nodes [20]. 


Similarly, an individual i that is susceptible at time t becomes infected at a rate !3f [ that 
depends on the state of the entire population X{t). This is known as an exogenous transition 
because it is influenced by factors external to the node itself. We discuss this at length in the 
sections to come. Figure [2] shows the simple interaction described above. 


Remark 1 (Other spreading models) We note here that this article excludes chain binomial 
models (e.g., Reed-Frost model) and other similar types of models from percolation theory. 
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Depending on the application at hand, the model for the spreading dynamics can vary. The 
main difference in the models we consider and models like Reed-Frost are that we want to 
allow infected individuals to continuously try to infect healthy ones. In the Reed-Frost model, 
an infected person only has one chance of infecting a healthy person. However, when thinking of 
a virus like the flu, a healthy person is continuously in danger of becoming sick when in contact 
with an infected individual, rather than a one-time chance. For instance, the Reed-Frost model 
might be more suitable for modeling the spreading of an email virus rather than an infectious 
biological disease. The interested reader should see [8], [211 for further details. 


We would also like to briefly point out that this other community is indeed active and 
working on similar types of problems as the ones we will highlight in this article. Many works 


exist on forecasting the cascading effects of a single infection or failure on a network [22], [231 


and how they can be mitigated through vaccination [24 ]. Conversely, one may be more interested 
in finding the most influential nodes or where to start an infection in a network to reach as many 


people as possible [25], [26]. This is often referred to as a seeding problem. Further extensions 


study attack and vaccination strategies on these models [27] and even cases in which there are 


multiple contagions on multiple networks [28]. 


Classical models 


Based on the above discussion, the dynamics of the SIR model is described by a 3^- 
dimensional Markov process. The exponential size of the state space makes this model very hard 
to analyze. One standard method to simplify the analysis is to consider the evolution of the total 
number of healthy and infected individuals rather than the state of each individual separately. 


This is commonly referred to as population dynamics [29[, [ 301. Furthermore, the recovery and 
infection rates are often assumed to be the same for all individuals; that is, 8i = 5 and /3f ff = /3 eff 
for all i. The standard population dynamics assumes a well-mixed population which means all 
individuals affect and are affected by all other individuals equally. Figure [3] shows the described 
interactions of this well-mixed population. 
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Stochastic population models 


The SIR population model is described as follows. Letting N z (t), N R (t ) e {0,1,..., N} 
be the number of infected and removed individuals at some time t, respectively, the number of 
susceptible individuals is necessarily given by N s (t ) = 1 — iV 7 (f) — N R (t). A common choice for 
the infection rate is given by /3 eS = /3A r/ A r ‘ s [7|, [311, [32[ for some f3 > 0, known as the mass 
action law. In other words, the rate at which the total number of susceptible individuals become 
infected is proportional to the product of the number of susceptible and infected individuals in 
the population. The state at some time t + At is then given by 

(N 1 + 1, N r ) with probability j3N I N s At + o(At), 

(N 1 — 1, N r — 1) with probability A/V 7 Af + o(Af), (1) 

(N s , N r ) with probability 1 — + <5 iV 7 ) At + o(At). 

For the SIS model we simply force N R = 0 at all times which simplifies this to 
N 1 + 1 with probability /3N I N s At + o(Af), 

N 1 — 1 with probability 5N 1 At + o(At), (2) 

N 1 with probability 1 — + 5 A 7 ) At + o(At). 


(N 1 ,N r ) 


N 


Removing the explicit definition of time, the SIS process can then be seen as a random 
walk on a line [32]-[35] (a similar Markov chain can be described for the SIR model) 


N 1 -P- N 1 + 1 with probability 
N 1 —^ N 1 — 1 with probability 


P{N- N 1 ) 
P(N — N 1 ) + 5’ 
5 

j3(N — N 1 ) + S' 


(3) 


An important observation about this model ([3]) is that it is a Markov chain with a single 
absorbing state N 1 = 0 in which all agents are healthy. In other words, once the entire population 
is healthy, the infection cannot suddenly reemerge. It is known from the theory of Markov chains 


that given enough time, the infection will eventually die out with probability 1 (see [361 for a 
review of Markov chains and relevant properties). Thus, the study of these systems is often 
interested in answering the question of when or how quickly the infection will die out. We 
comment on this later, in Remark [3} 
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To further simplify the problem, many works consider a deterministic approximation of 
these stochastic dynamics. In fact, the simpler deterministic dynamics we introduce next predate 
the introduction of the stochastic model above ©• 

Deterministic population models 


The models we present next are perhaps the two most studied epidemic models in the 
literature and is covered in a large number of different books |6], [8j, [10]-[f2], p7|-[|44|. These 


books also discuss a large variety of extensions including more complicated disease models that 
have more than two states, modeling for birth and mortality rates, different types of infection 
rates, and different categories for each disease state; for example, based on age or sex. We only 
present the most basic models here. 


Assuming a large population size N, we define p 1 = and p s = n - n ^- nR a s 

the fractions of infected and susceptible individuals, respectively. Then, we can write the 
deterministic SIR version of (jT]) as 


f = -/?//, (4) 

p 1 = /3p J p s — Sp 1 , 

and the deterministic SIS version of Q as 

p s = —^p 3 + dp 1 , (5) 

P 1 — /3p 7 p S — Sp 1 . 


These are derived by leveraging Kurtz’ theorem while assuming N to be very large [431. 


For simplicity, we will continue the analysis only for the SIS model but note that similar 
analysis can be done for the SIR model as well. Since the population size N is fixed and 
N R = 0, we have that p s = 1 — //, and the above equations ([5]) are redundant. Hence, they can 
be simplified to 


P 1 — /3p J ( 1 — p 1 ) — Sp 1 . 


( 6 ) 


Given an initial condition p / (0), this equation can be analytically solved [8], [45], [46|. The 
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solution is given by 


p\t) 


e (P~6)t 

P(e(P- S ') t -l) 

1=6 


/3t+ — 

n 1 


I- 

P 1 ( 0 ) 



P^s, 

p = 5. 


Given the exact solution of p [ (t), we are able to characterize its equilibrium points in the 
following result. 


Theorem 2 (Solutions to deterministic population model) The solution of p r (t) approaches 
1 — | as t —* oo for f3 > S, and 0 as t —» oo for (3 < 5. 


Remark 3 (Deterministic vs stochastic population models) It is important to note here that 
the deterministic models are only approximations of the stochastic models. A natural question 
to ask then is the following: What can the threshold result of Theorem [2] tell us about the 
stochastic model in ([3])? The first thing to note is that in the stochastic model, given enough 
time, the system will reach the disease-free state with probability 1. However, from Theorem [2] 
we see that for (3 > 5, the deterministic model will converge to an endemic equilibrium, meaning 
the disease never dies out. Thus, rather than studying the equilibrium values of the two models, 


the authors in [35], [45 [, [47] look at the expected time E[T] for the stochastic model to reach 
the disease-free equilibrium. Interestingly, they are able to show that for f3 < 5, the expected 
time E[T] is upper-bounded by -f 1 . On the other hand, when /3 > <5, the expected time E[T\ 
grows exponentially with N. • 


The analysis of the deterministic model results in a very precise threshold result that 
translates directly to the stochastic model as discussed in Remark [3] Threshold conditions are 
often given in terms of a reproduction number R 0 , which is the expected number of individuals 
a single infected individual will infect [7J, [48] over the course of its infection period. In other 
words, given a fully healthy population, if we infect person i at random, R 0 is the expected 
number of other individuals that will become infected over the course of agent i’s infection. 
This is a very useful metric with a critical value of Ro = 1. When R 0 < 1 the disease does not 
spread quickly enough, resulting in a decay in the number of infected individuals (in expectation). 


On the other hand, when R 0 > 1 the infected population grows over time (in expectation) [ 101. 
In the simple model considered above, the reproduction number is given by |. Furthermore, the 
exact solutions and asymptotic behavior of the system are easy to obtain. 
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The reproduction number is an important parameter that epidemiologists are interested in 
identifying for different diseases and environments [ |49| as it is a single number that can predict 
whether a certain outbreak of a disease will become an epidemic or die out on its own. Of 
course, the problem is that computing R 0 in general is not trivial, as there is no database for 
things like infection rates and recovery rates for different diseases. 


The main drawback of these population models is that they are very crude model derived 
by making many simplifying assumptions including (/) a homogeneous incidence rate (3 eff and 
recovery rate 5 for all individuals, (ii) a two-state model, (Hi) a constant population size, and (iv) 
a well-mixed population (or a contact network that is a complete graph). These drawbacks were 
very evident when scientists attempted to estimate the reproduction number of SARS in China 
in 2002-2003 and grossly overestimated it. This led to SARS scares making global headlines 
which eventually fizzled out because the actual reproduction number was far less than estimated 


due to the crude population models. More details on how this occurred can be found in [50], 
but the upshot is that more refined models are needed. 


Network models 


In order to create more refined epidemic models, it is clear that we cannot simply lump 
an entire population into two compartments defined by a single number. Ideally, we would be 
able to model the states of all N individuals independently and allow for arbitrary interactions 
among them. Not surprisingly, this is not a trivial task. 


In this section we are interested in spreading processes on a given, arbitrarily topology. 
Before jumping into the models we are interested in, we should mention that we are skipping 
over a plethora of work that has been dedicated to extending the population models to structured 
network models. More specifically, before jumping to completely arbitrary networks, there is a 
large body of work that studies various, specific structures. For instance, some works study 


how a disease spreads on a two-dimensional lattice or star graph [51 ]—[53]. Others consider 
more complex interconnection patterns, such as power-law and small-world networks, that still 
have some exploitable structure. In this context, a common method to analyze these networks 
is to assume that nodes are infected at a rate proportional to the number of neighbors they 


have [54[-[59[. These methods are justified depending on the assumptions enforced on the 
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network topology. A review of these types of models can be found in (60J. In what follows, we 
present epidemiological models on arbitrary networks. 


Stochastic network models 


Here we describe an epidemic model described as a continuous-time networked Markov 
model. Consider a network of N nodes represented by a connected, undirected graph Q = (V, E) 
where V is the set of nodes and E C V x V is the set of edges. We define the matrix A 6 n 
as the adjacency matrix of the graph, defined component-wise as a t] = 1 if node i can be directly 
affected by node j, and a V) = 0 otherwise. See “Graph Theory” for further details. 

We let Xi(t) denote the state of node i at time t, where Xft) = 1 indicates that i is 
infected and Xft) = 0 indicates that i is healthy at time t. Infected nodes can transmit the 
disease to its neighbors in the graph Q with rate 8 > 0. Simultaneously, infected nodes recover 
from the disease with rate 5 > 0. Figure [4] shows the described interactions on an arbitrary 
network. One can model this spreading process using the Markov process 

x i : 0 ->• 1 with rate 0 EjeVi X E (?) 

Xi : 1 —>• 0 with rate 5. 

Notice that there exists one absorbing state in this Markov process (corresponding to the disease- 
free equilibrium) that can be reached from any state X(t) = [X\ (t),..., X^(t ')\ r . This implies 
that, regardless of the initial condition X(0), the epidemic will eventually die out in finite time 
with probability 1. 

An interesting measure of the virality of a spreading process is the expected time E[T] 
it takes for the epidemic to die out. In |6TJ and [ (62| , we find the following threshold conditions 
in terms of the infection strength, defined as r = |. 


Theorem 4 (Threshold for sublinear expected time to extinction) If 


T < 


\nax(A) 

where A max (A) is the maximum real eigenvalue of A, then 

log N T 1 


E[T\ < 


$ ~ /3A max (A) 


for any initial condition X (0). 
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Note that Theorem [4] only provides a sufficient condition for ‘fast’ extinction of a disease. 
There have been many efforts to determine whether this condition is also necessary, but at the 


time of writing this remains an open question on general graphs. The works [43], 
show that there exists some critical value r c of the infection strength for which the expected 
time to extinction grows exponentially with N when r > r c . The following result formali z es this 


statement and provides a lower bound on the critical values [65], [66 [; however, we should also 
note that stronger statements exist when considering graphs with a fixed structure (e.g., lattice, 


star) [431. 


Theorem 5 (Threshold for exponential expected time to extinction) There exists 

1 

T r > 


^max {A) 

such that for r > r c , the expected time to extinction E[T] = O (e kN ), where k depends on r 
and the structure of the graph Q. 


The maximum eigenvalue A max (Al) of an adjacency matrix is a parameter that captures 
how ‘tightly connected’ the graph is. In general, more connections means a large A max (Al). 
Intuitively, the results of Theorems [4] and [5] are saying that the more tightly connected the graph 
is, the easier it is for a disease to spread. 


It is worth mentioning that although the result of Theorem [4] provides an upper bound 
on the expectation of the extinction time, the possibility of a persisting epidemic is not ruled 
out. For example, it has been shown for star graphs that regardless of the infection strength r, 
there is a positive probability that the time to extinction is super-polynomial in the number of 


nodes [611, [67], [68]. Furthermore, for high-degree or scale-free networks (such as preferential 
attachment [68] or power-law configuration model graphs [|61]), it has been shown that this 
threshold goes to zero as the number of nodes increases j69[ because the maximum eigenvalue 
grows unbounded with N. 


Deterministic network models 


Here we describe the deterministic model of the SIS dynamics over arbitrary net¬ 
works [ j20| , [70]—[74 ]. We begin by assuming homogeneous recovery and infection rates, although 
this assumption will be relaxed in the following section. The natural recovery rate of each node 
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is given by S > 0 and the infection rate at which a node is affected by infected neighboring 
nodes is fd > 0. The dynamics of the spread is described by the set of ordinary differential 
equations 

N 

pi = -Spi + ^2 aijPpjtf. - Pi ), (8) 

3 = 1 

where p t (t) G [0,1] describes the (approximated) probability that an individual i is infected 
at time t. See “Networked Mean-Field Approximations” for further details. This variable has 
another interesting interpretation in the context of metapopulation models. In a metapopulation 
model, each node does not represent an individual, but a large subpopulation (such as an entire 
district or city). In this context, pi can be interpreted as the fraction of the i-th subpopulation 
that is infected. See “Meta-Population Models” for further details. 


As with all other epidemic models, we see that the disease-free equilibrium Pi = 0 for 
all i G A’} is a trivial equilibrium of the dynamics. We are now interested in finding 

conditions such that this equilibrium is globally asymptotically stable. Letting p — (pi,... ,Pn) T 
and recalling the infection strength r = f, the following result from [74]-[77] characterizes the 
convergence properties of these dynamics. 


Theorem 6 (Threshold condition for networks) Given the dynamics ([8]) for any p(0) 0, the 

equilibrium p* = 0 is globally asymptotically stable if and only if t < t— W-. Furthermore, for 

^max\A) 

t > A there exists p** such that p** is globally asymptotically stable. 

Remark 7 (Deterministic vs stochastic network models) Similar to our discussion in Re¬ 
mark [3j there is a connection between the deterministic result in Theorem [6] and the stochastic 
result in Theorem [4] Since X = 0 is an absorbing state, the stochastic dynamics will eventually 
reach the disease-free state with probability 1. However, Theorem [6] claims that for /3A max (A) > 5 
the deterministic model will converge to an endemic equilibrium, meaning the disease never 
dies out. To resolve this apparent contradiction, we again focus our attention on the expected 
time E[T] for the stochastic model to reach the disease-free equilibrium. Remarkably, Theorem [4] 
provides a sufficient condition for a disease to quickly die out that is in agreement with the 
threshold result of Theorem [6} However, as suggested by Theorem [5j it has not yet been shown 
whether the same threshold condition holds for persistence of the disease in the stochastic model. 
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A major drawback of the above model is that it assumes a constant infection rate f3 and 
recovery rate 5 for all individuals. To further refine the model, we are interested in allowing 
different recovery rates for each person and different infection rates for each type of contact. This 
allows for a much more general model that can capture more realistic scenarios. For instance, it 
is not fair to assume that everyone you come in contact with has an equal chance to infect you. 
A family member or a spouse is much more likely to infect you than a casual acquaintance. To 
capture these heterogeneous effects, we develop heterogeneous network models next. 

Heterogeneous network models 


Here we describe the dynamics of the SIS model with heterogeneous recovery and 
infection rates over arbitrary strongly connected directed graphs Q = (Vi E). The recovery 
rate of node i is given by Si > 0. In our exposition, we consider an edge-dependent infection 
rate. In other words, the infection rate at which a node i is affected by an infected node j is 
given by > 0 if (i,j) G E. For simplicity, we let fj tJ = 0 if (?', j) £ E. The dynamics of the 


SIS model in an arbitrary network is described by [751 


N 


Pi = -SiPi + ^2/3 ij p j { 1 - Pi), 

3 = 1 


(9) 


where p, £ [0,1] can be seen as either the fraction of the i-th subpopulation that is infected (in 


the metapopulation case), or the probability that an individual i is infected [73], [75]—[79[. 


In this model, the disease-free equilibrium p t = 0 for all i e {1,..., A r } is again a 
trivial equilibrium. In what follows, we derive conditions for this equilibrium to be globally 
asymptotically stable. Let p = (pi,,P m)' denote the state vector of the system, D = 
diag(<5i,... ,5n) the diagonal matrix of recovery rates, and B = [/3ij\ the matrix of infection 
rates. The dynamics ([9]) can then be written as 


p = (B - D)p + h, 


where hi = fiijPiPj- The following result from [75], [76], [80] characterizes the 

convergence properties of these dynamics. 


Theorem 8 (Threshold condition for heterogeneous networks) Given the dynamics in 
for any p(0) f 0, the equilibrium p* = 0 is globally asymptotically stable if and only if 
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A max{B — D) <0. Furthermore, for A max (B — D) >0, there exists p** e »S,i) such that 
p** is globally asymptotically stable. 


These stability results have been recently extended to a number of more complicated 
models such as the three-state SAIS model (81J, the four-state G-SEIV model []82], and even 


the SI*V* model with an arbitrary number of states [83]. 


Now that we have provided a basic understanding of how the SIS process evolves and the 
connections between the stochastic processes and their deterministic approximations as discussed 
in Remark |7] we are now ready to formulate and study some relevant control problems. 


Control of Epidemics 


In the previous section we presented several approaches for modeling the dynamics of 
spreading processes taking place on arbitrary contact networks. We have also analyzed these 
models and introduced several stability results for both the deterministic and stochastic cases. 
In this section, we describe several results aimed at controlling the dynamics of the spreading 
processes. 

Ideally, we are interested in controlling the stochastic network models to stop the 
spreading of a disease as quickly as possible. However, before getting to the details, we must 
begin by talking about our effective ‘control levers’ in treating an epidemic. For simplicity, let 
us consider the heterogeneous SIS dynamics ([9]) 

M 

Pi = -SiPi + ^2 PijPji 1 - Pi), 

3 =1 

as a meta-population model with M subpopulations. That is, each node i is some subpopulation 
(such as a town) of n % individuals in a larger population (such as a country) of N individuals 
(see “Meta-Population Models” for further details). The parameters we have to play with are 
then the recovery rates Si for each subpopulation and the infection rates fiij that describe the 
interactions between various subpopulations. 

In order to mitigate the effects of an epidemic, in general we would like to increase the 
recovery rates Si and decrease the infection rates (i l3 . Increasing the recovery rate of a given 
subpopulation can be done by providing better treatment to sick individuals. For instance, by 
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allocating more resources to this particular subpopulation they can afford more doctors or better 
methods of treatment for fighting a particular disease. Decreasing infection rates can be done in 
numerous ways. Limiting traffic/travel between subpopulations can help decrease the infection 
rate. Completely quarantining a subpopulation i is equivalent to setting (3ji = 0 for all j since i 
can no longer affect other subpopulations. Other ways of decreasing infection rates include 
milder methods of prevention, such as distributing masks to a population to minimize chance 
of infection; or even simply raising awareness about a disease to make people less likely to 
contract the disease. 

Clearly if we had infinite resources, treatment power, and simply quarantined everyone 
the disease would likely die quickly; however, this is not a feasible solution. Thus, given a 
fixed budget of some sort, it is imperative to identify which parameters specifically are the most 
important in order to mitigate the effects of the disease as much as possible. We formulate these 
problems and discuss the current state of the art next. 

Spectral control and optimization 

Here we are interested in various optimal resource allocation problems. More specifically, 
given a fixed budget, the idea is to optimally invest resources to best hinder the spreading of 
a disease. Leveraging the results of Theorems |4jj6} a natural option to mitigate the effects of a 
possible epidemic is to make A max (B — D ) as small as possible. 

We first discuss the homogeneous SIS dynamics ([8]) where 5 and 3 are fixed parameters. 
Hence, we are interested in making A max (/1) as small as possible. This can be achieved by 
modifying the network structure. The work [84] studies the effects that the network structure 
has on this maximum eigenvalue. In this work, the authors study how to decrease A max (A) in 
one of two ways. The first is to remove nodes from A. This might physically be done by either 
quarantining or immunizing certain individuals, making them unable to contract the disease and, 
more importantly, to spread it. Another way to reduce A max (/l) is to remove links rather than 
completely removing nodes. This might physically be done by limiting traffic between certain 
cities or limiting interactions between certain individuals. The caveat is that we are interested 
in doing this while removing the least amount of nodes or edges since these actions are likely 
quite costly in the real world. The node and link removal problems of interest are then described 
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as follows. 


Problem 9 (Optimal node removal) Given an original graph A and a fixed budget C' > 0, 
minimize X m ax(A) by removing at most C nodes from A. 

Problem 10 (Optimal link removal) Given an original graph A and a fixed budget C > 0, 
minimize \ max (A) by removing at most C links from A. 


Unfortunately, the node and link removal problems described above are NP-complete 


and NP-hard, respectively [85]. As a result, several papers instead solve convex relations or 


propose heuristics to approximately solve these problems. A simple example is one in which 
the nodes with the highest degrees (largest numbers of neighbors) are removed one by one 
until the budget is exhausted. Other heuristics are based on various network metrics, such as 


betweenness centrality [86], PageRank [87], or susceptible size [88], to decide which nodes 
should be removed first. Similarly, there are works that are concerned with link removal rather 
than node removal [[85|, [89], [9(1. In [|9l], the authors solve a convex relaxation of the problem 


and effectively project its optimal solution onto the original problem. 


Unfortunately, the authors in [92| study the worst-case scenarios of these suboptimal 


strategies to show that network-based heuristics can perform arbitrarily poorly. Thus, it is hard to 
evaluate a priori how well a suboptimal solution to Problems [9] and p~Q] will perform. Furthermore, 
completely removing nodes or even links might not be feasible solutions anyway as this would 
require fully quarantining certain subpopulations or completely shutting down certain roads or 
methods of travel between various subpopulations. 

Instead, let us now consider the heterogeneous network model in ([9]) and tune the values 
of the parameters <5* and ffij rather than completely changing the network structure. The authors 


in [931 formulate this as an optimization problem to minimize the steady-state infection values 


over heterogeneous recovery rates. A gradient descent algorithm is then proposed to find feasible 
local minima solutions. Another alternative is to utilize the result of Theorem [8j In this direction, 
several papers consider the minimization of A max (7i — D) under various constraints. The effect 
of minimizing this eigenvalue is to maximize the exponential decay rate of the system towards 
the disease-free equilibrium. 


While tuning the spreading and recovery rates, one can consider a discrete optimization 
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setup in which one can only tune these rates within a discrete set of feasible values. This problem 
has been shown to be NP-complete in (94) . Alternatively, one can consider a relaxation in which 
these rates can take values in a feasible continuous interval. In this case, the works in 


196] propose efficient methods for allocating resources to minimize the dominant eigenvalue 


of relevant matrices. In [97] and [98], the problem of minimizing A max (T> — D) is cast into a 


semidefinite program framework for undirected networks. In [99[, [100], this problem is solved 
for directed graphs using geometric programming where the solution can be obtained using 
standard off-the-shelf convex optimization software. Furthermore, geometric programs allow for 
the simultaneous optimization over both the infection rates and recovery rates. See “Geometric 
Programming” for further details. 


In what follows, we present a simplified version of the optimization problem considered 


in [100] and show how it can be efficiently reformulated as a geometric program. Consider the 
deterministic heterogeneous SIS model ([9]) with natural recovery rates Si = !u > 0 and infection 
rates d, = 3 r > 0 for all i G {1,..., N}, where 3ij = 3i for j G Af ' n and j3 tJ = 0 otherwise. In 
other words, the rate at which a node i is infected is a node-dependent parameter rather than an 
edge-dependent one. We then assume we are able to pay some cost to increase S z up to some 
maximum Si > Si- Alternatively, we can also pay a cost to decrease 3i down to some minimum 
0. < 3 t . The control parameters are then given by 6, and 0 h where < 5, < 3 ■ < 3i 3 3 r - 


Assume we have access to cost functions fi(Si) and gfiPp describing the associated cost 
to increase Si and decrease 3 t , respectively. In this context, given a fixed budget C > 0, our goal 
is to minimize A max (T> — D) while satisfying the constraint that the total cost does not exceed 
the given budget. This problem is formally stated below. 


Problem 11 (Budget-constrained allocation) Given a fixed budget C > 0, 

minimize A max (-E> — D) 

(A A )£i 

such that Y3i=\ fi(Pi) + 9i(Si) < C 
Pi < Pi < Pi 
Si<Si< Si. 

Note that solving Problem |TT] is not trivial since the objective (maximum eigenvalue) 
function is not convex in general. However, the following result guarantees that, under mild 
assumptions on the cost functions, this problem can be solved exactly by rewriting it as 
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a geometric program which can be solved using standard off-the-shelf convex optimization 
software. See [100] for further details on this equivalence. 


Theorem 12 (Solution to budget-constrained allocation problem) Problem m can be solved 
by solving the following auxiliary geometric program 

minimize A 

r —' i n 

5 't'U'i j ^i 

such that Y^j =1 a ijfii u j + $i u i < A Ui, 

E UWj)+m)< c , (10) 

(ft -Si < Si < (ft - Si, 

a < a < 

for all i G (1,..., N}, with (ft > rnaxj Sj and gj(Sj ) = gj((ft ~ $j)> where (3* and 5* = (ft — 6* 
solve Problem [77] with rate A max (Q) < A* — 0. 


We have only considered node-dependent infection rates $ here. This is extended to the 


case where these rates can be controlled over edges in [101]. Aside from the discussed SIS 
model, other recent works have also applied these ideas to more general models. The authors 


in [1021 formulate the SDP problem for a three-state SAIS model proposed in [81] in which 


alertness to a possible epidemic is also modeled. A general four-state SEIV model is considered 


in [1031 for which the authors develop equivalent geometric programs to optimize the dominant 
eigenvalue over various parameters of the model simultaneously. 

These types of optimal allocation strategies have been recently compared to fair strategies 


in [ 104], where resources must be allocated evenly across all nodes, to show their effectiveness in 
targeting resources rather than evenly spreading them. However, there are still some drawbacks 
of these spectral control approaches that need to be addressed before we can fully take advantage 
of their solutions in weakening the impact of diseases in the future. 


The first main drawback of these approaches is that they do not take into account the 
current state of the system. This means that even nodes that are not at immediate risk of being 
infected might be allocated resources to raise their recovery rates or decrease their infection rates. 
Second, solving these problems exactly requires a lot of knowledge. In addition to knowing the 
natural recovery rates and infection rates, we have assumed that we have exact knowledge of the 
entire graph which is a bit of a stretch. Third, these are centralized solutions which may take 
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a long time to compute. Although we have been able to solve some variants of this problem 
as discussed above efficiently (in polynomial time), this still may not be fast enough if these 
networks are very large. Lastly, we have also assumed that once the optimal solution is found, 
we are able to instantaneously set the recovery rates and infection rates to the desired values. 

We will discuss the current efforts on how each of these problems are currently being 
addressed and what still needs to be done in the sections to come. We begin in the next section 
by relaxing this first issue by looking at optimal control problems with feedback, rather than 
one-time optimal resource allocation solutions. 

Optimal control 

Here we discuss various optimal control problems formulated for mitigating epidemics 
under the SIS and SIR dynamics. However, since very little work has been done for the network 
models thus far, we start by looking at the classical models. 

Classical models 


We begin by recalling and slightly modifying the SIS population model ([6]) to account 
for a control action. Following ]105 1, we rewrite the original SIS population model with 8 — 


P = ftp (1 -p ) - 8 x p , 


( 11 ) 


where 8i > 0 is the natural recovery rate of an individual. We now assume that we are able to 
control this system by increasing the recovery rate of individuals in the population from Si to 
82 > 81 . This can be achieved, for instance, by allocating antidotes or providing other forms of 
treatment to a fraction of the population. Our control signal u £ [0,1] is then the fraction of the 
population that we provide treatment to. For simplicity, we assume that we are able to instantly 
affect the recovery rates of any number of individuals in the population. The dynamics of the 
controlled SIS population model is then given by 

p 1 = Pp\ 1 - P 1 ) ~ ((1 - u) 8 i + uS 2 ) p 1 ■ (12) 

Applying the result of Theorem [2j we obtain the following corollary for a fixed u{t) = a. 
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Corollary 13 (Population dynamics threshold condition) The solution of p 1 (t) approaches 0 
as t —> oo for 

P-5 1 

u - a - T- 

Since u G [0,1], Corollary [73] implies that if (j3 — 8f)/(8 2 — 5i) > 1, the disease is too strong 
and will never die out regardless of our control choices. On the other hand, when 8\ > /3 the 
natural recovery rate is high enough to ensure extinction of the disease without any control action 
(u = 0). Otherwise, it should be easy to see now that if we want to exert a minimal amount 
of a fixed control signal to ensure extinction of the disease, we can simply choose u = 3 / . 
However, it may be the case that we are still willing to use control such that the infection dies 
out faster than it would naturally. For instance, having a population with many sick individuals 
could incur a drastic social cost that could be instead offset by a smaller cost of treatment. We 
formulate this tradeoff as an optimal control problem next. 


Let the cost of treatment be linear with the number of individuals treated, and similarly, 
let the cost of infection be linear with the number of infected individuals. We are then interested 
in minimizing the objective function 


Jt — (cp 1 (f) + du(t))dt, 


(13) 


where c > 0 is associated with the cost of infection, d > 0 is the associated with the cost of 
treatment, and T > 0 is the time horizon. Utilizing Pontryagin’s maximum principle, it can be 


shown 105 ]—[ 107) that the optimal solution is given by 

f {0} for f(t ) > 0, 

u *(t) e l [0,1] for/(f) = 0, 

[ {1} otherwise, 

with 


/(f) = i’P 1 ^ - <5i) + d, 
where is the costate variable with dynamics 

ip — c ip(P( 1 - 2 p 1 ) - ((1 - u)5i + u8 2 )). 


It can now be shown [105] that for [ 5/(82 — 3i) < c/d, the optimal solution is to initially 
treat the entire population until some time t' at which nobody should be treated. For [ 8 /( 82—81 ) > 
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c/d, the optimal solution is u(t) = 0 for all t G [0, T], This bang-bang solution with at most one 
switch is very common in similar problems. Other works with this same kind of solution have 


been studied in many different variations of this problem that consider efficiency of control 1081 


or control over both 5 and f3 simultaneously 1091. Other models have also been considered 
such as the SIR model [ 11 Oj ] with different incidence rates [111], [l ~j~2] | or a four-state SIRD 


model [ 1131. 


Although the bang-bang solution is common, it is possible to obtain other types of 
solutions for different formulations of the optimal control problem. For example, it is shown 


in 114], [ 1151 that for alternative problem formulations, the optimal solution may not be a 


bang-bang controller for certain classes of cost functions. In [116], an SIR model with quadratic 
control costs over both 5 and (3 is considered. In this case, the optimal solution is again not 
a bang-bang controller. A four-state SIRC model for which the optimal solution is again not a 


bang-bang controller is considered in [117]. 


In the problem above, we have assumed that we are able to change the control signal 
instantaneously. Other works consider the case in which the rate of the control signal (its time 
derivative) can be controlled instead [33], [34], [|TT8j], [119]. We omit the technical details of 
these works as the methods are very similar to the example presented above. Interestingly, the 
results from these works often admit bang-bang controllers with at most one switch as optimal 
solutions as well. As a final note we acknowledge that in different contexts one may be interested 
in the problem of maximizing the impact of a spreading process (for instance a viral marketing 


campaign) [ 1121, 1201 rather than minimizing it. 


Network models 


As mentioned before, the population models are quite crude in general as they lump 
an entire population’s state into just a few numbers. Thus, we now turn to optimal control of 
networked models but note that, so far, very little work has been done in this realm. Three 


relevant papers that consider this problem in the context of networks are [107], [ 1211, [1221. We 
start our exposition by proposing an optimal control problem for SIS dynamics that has yet to 
be solved. 
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Recall the SIS network dynamics with heterogeneous recovery and infection rates <(9]), 

N 

Pi = + Pi)- (14) 

j =i 

From Theorem [8] we know that a necessary and sufficient condition for extinction is A max (7i — 
D) < 0. In the previous section, we utilized this result as a constraint to solve optimal allocation 
problems. Instead, we are now interested in solving optimal control problems, where rather than 
solving a one-time optimization to determine the curing rates 5 t , we allow them to vary over 
time. 

Problem 14 (Optimal control of an SIS network) Given a linear cost of infection Ci and 
control di for all i G (1,..., N}, minimize 

'Y] CiPi(t) + diSi(t) \ dt (15) 

subject to the dynamics Cl and Si(t ) G [<5, 5] for some 0 < 6 < S for all t G [0, T\. 



This, along with most of its variations, is currently an open problem. Variations include 
problems similar to the optimal control problems for deterministic population models discussed 
earlier: control over infection rates, non-instantaneous control, different objective functions, etc. 


The only work we are aware of that has tackled this problem is 1221, where the authors study 


the linearization of ( fT4] ) around the disease-free equilibrium. The authors are able to show for 
the linear dynamics that the optimal solution is a bang-bang controller with at most one switch, 
similar to many results obtained for the population models. 


Although Problem 14 is still an open problem for the SIS dynamics, a closely related 


problem has been successfully solved in the context of containing computer viruses [107], 
[ |121| . Here we present a simpler version of the problem originally posed in [ 1071. Consider the 
dynamics 


Pi = 


N 

-p s i 3i ’i\ 


j - pfpSi, 


3 = 1 


N 


Pi = Pi ~ P i niR ‘ 

3 = 1 


Ui 


(16) 


pf = PiPfui + ptmRiUi, 
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where as before pf and pj are the fraction of a subpopulation that are susceptible and infected, 
respectively. Then, pf = 1 —pf — pi is the fraction of individuals that are ‘removed’. This refers 
to individuals that are immune from the infection, whether this means they were vaccinated or 
recovered from the disease and are no longer susceptible to it. Additionally, Ui is the control 
that dictates the rate at which susceptible and infected individuals become removed. 


Problem 15 (Optimal control of an SIR network for malware epidemics) Given a linear 
cost of infection Ci, control h] and lif, and benefit of recovery l t for all i G (1,..., A’}, we 
would like to minimize 


N 


Jt = 


Y + Ci Pi + P? h l U i + P?(Pi + Pi) h i u i dt i 


(17) 


, 2=1 


subject to the dynamics © and lift) G [0, for some Ui > 0 for all t G [0, T\. 


The following result follows from Pontryagin’s maximum principle [107]. 


Theorem 16 (Optimal control of an SIR network for malware epidemics) There exists 
Ti G [0, T] for all i such that the optimal control is given by 


u*(t) 


Ui for t < Ti, 
0 for Ti <t <T. 


Again, Theorem |T6] is consistent with many other optimal control solutions for epidemics 
in that the optimal solution is a bang-bang controller with at most one switch. Given the 
recentness of these results, there are still lots of variations of this work that need to be studied. 
Although the dynamics ( fl6] ) considered here is very similar to the epidemic models we have 
discussed throughout the article, it is not immediately applicable due to the term R,u t . In the 
context of patching, R is a state of nodes who have a ‘patch’ and are thus immune, and so they 
can spread this patch to healthy and infected nodes. However, this concept does not translate 
directly to general epidemics; a sick person cannot get better by interacting with healthy people. 
We expect these types of problems to be solved for epidemics in the very near future. 

In many of the problems we have discussed above, we have assumed that we have direct 
control of the infection rate fiij and the recovery rates d t . However, this simplistic scenario 
assumes that we can control these parameters for the entire population instantaneously which 
is unfeasible in the context of disease spreading. In an effort to address this oversimplification, 
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there is a rising body of current work in which more realistic ‘control’ actions are explored. We 
discuss these next. 

Heuristic feedback policies 


Here we briefly review various models that are used to capture possible human behaviors 
or other countermeasures employed to deter the spreading of a disease. Rather than explicitly 
attempting to control the SIS dynamics as described above, the works we discuss here are 
essentially extensions to the SIS model for which stability conditions are derived. The models 
are created by assuming various actions people might take, and then the ‘closed-loop’ system 
stability is analyzed. More specifically, rather than separately considering a model and control 
strategies, the model and control strategies are co-developed to yield a sense of ‘closed-loop’ 
control model. For lack of better terminology, we refer to these as heuristic feedback policies. 


More specifically, many works consider various feedback strategies that determine when 


nodes or links should temporarily be removed 891, 123]—[ 129]. Closed-loop models are then 
constructed for the various strategies whose stability properties can then be analyzed. These 
strategies are generally based on some sort of perceived risk that individuals have of becoming 
infected, causing them to either remove links to infected neighbors or completely remove 
themselves from the network (for example, by staying home from work or becoming vaccinated). 
We begin by looking at the simpler classical models, then later show how these can be extended 
to network models. 


Classical models 


As mentioned above, these so called heuristic feedback policy solutions are all essentially 
different epidemic models for which stability results are obtained. As an illustrative example, 
we consider the work ]128 1, where in addition to the susceptible state S and infected state I, 
an additional protected state P is introduced. The protected state refers to individuals who have 
decided to immunize themselves in one way or another, and are thus not immediately susceptible 
to contracting the disease. The model is described as follows. Letting Y, be the number of 
infected in-neighbors a susceptible node i has in an appropriate graph, node i transitions from 
the susceptible state S to the infected state / with rate 8 Y- However, a node in the protected 
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state P transitions to the infected state with rate l3 0 Y t where (3 0 < f3 captures the decreased 
risk of infection due to being ‘protective’ or ‘alert’. A type of control is then to decide how 
susceptible individuals transition to the protected state. Finally, as in the normal SIS model, 
individuals that are infected naturally recover to the susceptible state with a natural recovery 
rate 5. Figure [5] shows the interactions of this three-state SPIS model. The authors then consider 
the extension of the SIS population dynamics ([5]) (by assuming a complete network topology 
meaning all individuals are equally likely to affect one another) to include the protected state 
given by 


p s — —^p 3 + dp 1 — p s f(p s ,p ! ,p p ) + p p g{p s iP 1 ,p p ), 

p 1 — Pp r p S — dp 1 , (18) 

•P Q n/ Q T P\ p / e T p x 

p = p j(p ,p ,p ) - p a{p ,p ,p ), 


where /(■) and g(-) are functions that determine how susceptible individuals are protecting 
themselves. Recall that p s corresponds to the fraction of individuals in a population that are 
in the susceptible state with p 1 and p p defined similarly for the infected and protected states, 
respectively. As in the case of the deterministic SIS population dynamics ([5]), one of these 
equations is redundant and can be removed by using p s + p 1 + p p = 1 because we assume 
a constant population size. We refer to this as the three-state Susceptible-Protected-Infected- 
Susceptible (SPIS) model. 

The authors then explore different strategies for designing / and g, and analyze the 
stability of the system for these choices. As mentioned above, we refer to this as a ‘heuristic 
feedback policy’ because a specific control structure is already defined and built into the model, 
rather than the objective of the work to be designing the controller itself. More specifically, in 
the example above, if we are free to choose the functions /(•) and g(-) arbitrarily, it is clear that 
the best thing to do is simply set g(-) = 0 and have /(■) be very large. This means everybody 
simply protects themselves very quickly, for which it is easy to imagine that the disease will die 
out quickly as well. Instead, it is useful to explicitly model a cost for infection and/or control 
as we did in the previous section. Next, we show how these models can be extended to network 
models. 
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Network models 


Similar to the classical models, here we discuss models that are used to capture possible 
human behaviors or other countermeasures employed to deter the spreading of a disease on 


networks [130], [[1311- As before, many works consider various feedback strategies that determine 


when nodes or links should temporarily be removed [20], [ 1321—[ 135]. Closed-loop models 
are then constructed for the various strategies whose stability properties can then be analyzed. 
These strategies are generally based on some sort of perceived risk that individuals have of 
becoming infected, causing them to either remove links to infected neighbors or completely 
remove themselves from the network (for example, by staying home or becoming vaccinated). 


As an example, consider again the three-state SPIS model ( p~8j ) presented in [ 1281. 
However, we are now interested in the network version of the population dynamics there. To do 


this, we create a model very similar to the three-state SAIS model presented in [81] where the 
authors introduce an alert state A which is similar to the protected state P we consider. This 
state captures the possibility of human behaviors and actions lowering the chance of contracting 
a disease. For simplicity, we consider homogeneous parameters so the recovery and infection 
rates are set the same at all nodes. The deterministic version of this model is then given by 

N 


Pi = -Ppf a v pI ] + 5 pI - pffi(p S :P T ,p P ), 


3 = 1 


N 


N 


pi = Ppf a vPj + a >.Pi ^Pj - 6 Pi’ 


3 =1 


3 = 1 


N 


Pi = Pi fi{p S ,P ! ,P P ) - PoPi a ijPj ’ 

3 = 1 

where /j(p 5 ,p 7 ,p p ) is a function that determines how susceptible individuals are protecting 
themselves. Conditions can then be derived for the parameters and the function / such that the 
disease-free equilibrium is globally asymptotically stable [811. The authors in [ 136] then treat 
the design of this function / as an optimal information dissemination problem. However, as in 
the population dynamics case, these are very structured methods of ‘control’ that ultimately get 
built into the models. 

A very large shortcoming of these types of solutions is that they are too specific. A 
very specific model with a specific control structure is proposed and studied. Unfortunately, it is 
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often unclear what type of spreading process each model is good for describing, if any. In [137], 
Hethcote does a great job highlighting this fact that there are far too many slight variations of 


existing models. We close this section with a small anecdote from the epilogue of [137], after 
effectively proposing 10 5 different models, to make us really think which models are actually 
interesting: 


“In the book [sic], A Thousand and One Nights, Scheherezade had to entertain King 
Shahriyar with a new story each evening in order to avoid being killed. If they were 
mathematical biologists and she had only to present one new epidemiological model 
each night to entertain him, then she could have survived each night for at least 270 
years. Of course, the King would probably have become disenchanted by the “new” 
models if they were only very slight variations on previous models and would have 
killed Scheherezade. Similarly, referees (the Kings) might become disenchanted if the 
papers which they receive contain models which are only slight variations on previous 
models. Thus I suggest that we as modelers and mathematicians should be cautious and 
not assume that every mathematical analysis of a slightly different model is interesting.” 

We review this issue and other technical challenges in the following section. 

Future Outlook 


In the previous section we have provided a high level overview of the current state of 
the art involving the control of epidemics. However, there are still substantial shortcomings of 
the results presented that need to be taken account for to take full advantage of their proposed 
solutions. Here we highlight several of the main research challenges, how they are currently 
being addressed, and what still needs to be done. 

1) All control methods we have discussed so far have been for deterministic models. 

Ideally, we are interested in ultimately controlling the original stochastic epidemic models from 
which the deterministic models are derived. Results like Theorems [4] and [5] help draw connections 
between the two, but these have only been done for simple cases so far. Furthermore, while these 
results help justify how using spectral control and optimization methods for deterministic models 
translate to the original stochastic models, it is unclear how the optimal control solutions found 
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for deterministic models relate to the stochastic ones. 


We are only aware of a few works that attempt to control epidemic processes on networks for the 


original stochastic models. In [671, the authors consider the SIS model with a simple heuristic 
control law where the curing rates Si for each node are proportional to the number of neighbors. 
The authors are then able to show that on any graph with bounded degree, this policy can 
achieve sublinear expected time to extinction with a budget proportional to the number of nodes 
in the network. A drawback of the above method is that it does not take into account the current 
state of the system. This means that nodes with no infected neighbors may be assigned high 


curing rates. Instead, the authors in [138] use a heuristic PageRank algorithm to allocate curing 


resources given a fixed budget, based on the initial condition of infected nodes. The authors are 
then able to provide probabilistic upper bounds on the expected extinction time. More recently, 
the work |139fl proposes an algorithm for which the expected time to extinction is sublinear 
using only a sublinear budget (in the number of nodes) for graphs satisfying certain technical 


conditions. In [1401, the authors consider a similar problem for which various algorithms are 
developed using a Markov decision process framework. 

2) All control methods we have discussed so far have admitted centralized solutions. 

This is a big problem since human contact networks can be massive in practice and it may not 
be computationally possible to solve these problems in a centralized setting. In this direction, 
distributed allocation and control strategies are an interesting alternative. Again, there are only 


a few recent works that have looked at this problem [107], [141 [-[143]. As more work in 
optimization and control of epidemic processes is being done, we desire and expect distributed 
versions of these algorithms to follow. 

3) All problems and solutions discussed so far assume no uncertainties. 

This is a very big issue in the context of epidemics. Throughout all the modeling, analysis, 
and control solutions we have presented so far, we have assumed perfect knowledge of 
everything including recovery rates, state information, and network structures. These are clear 
oversimplifications since in practice we would be more than lucky to have any of these parameters 
simply handed to us. A review of analysis and approximation techniques considering uncertainties 


in the spreading parameters is provided in [144], which is a large field of study in epidemiology. 
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In the context of control, far less work has been done for the case where the topology is 


unknown. In [145], the authors use observed infection data of a discrete time SIS process to 
estimate the network topology. Optimization and control methods can then be applied to the 
estimated topology; unfortunately, it is unclear how well these solutions will perform on the 
actual topology. Instead, a data-driven approach to optimally allocate resources has been recently 
proposed in [ 146[ , where only empirical data about the spreading of a disease is available. In 
this work, the authors assume that the spreading and recovery rates are unknown. Alternatively, 
the authors assume that the responsible health agency has access to historical data describing 
the evolution of the disease in a network during a relatively short period of time. In this context, 


the authors in 1461 propose a robust optimization framework to allocate resources based on 
historical data. 


Another large issue related to assuming perfect knowledge is assuming that we are able to 
set recovery and infection rates to whatever desired values we like. This is again a clear 
oversimplification and studies are needed into how various control solutions perform when these 
rates cannot be set exactly. 

The assumption of being able to observe exact state data is another big issue that has received 
little to no attention in the context of controls. These do not apply to the spectral optimization or 
heuristic feedback control methods, but are certainly important for the optimal control methods. 

4) We require much more general epidemic models. 


Although there has been a lot of work on modeling in addition to the SIS and SIR dynamics we 
have mainly focused on throughout this article, there is still a lack of generalized models. More 
specifically, a majority of works that study spreading processes begin with a single model with 
a fixed number of states and interactions. Many of these models are created by first looking at 


empirical data of a spreading process like AIDS [147] or a computer virus [69], then determining 
what type of model and how many states should be used to capture its behavior. Instead, few 
works propose much more general models with arbitrary numbers of states or layers on which 
the disease can spread 


|, [83], [148]. The further development, analysis, and control of these 
generalized models can allow rapid prototyping of models for spreading processes that might not 
even exist today, in addition to completely generalizing the myriad of specific models available 
today. 
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All models we have discussed in this article so far only consider the spreading of a single disease 
or process. Extending existing models to capture multiple diseases that co-evolve in a network 


has recently been gaining attention 149]—[ 1551. In general, these diseases are assumed to be 
mutually exclusive, meaning an individual can only be infected with one type of infection at a 
time. While it is discussed here in the context of disease and epidemics, these models are more 
aptly used in studying belief propagation or product adoption. For instance the mutual exclusion 
of infections is very relevant in competition in politics, such as Democrats vs Republicans, 
or competition in a marketplace, such as iPhone vs Droid vs Galaxy. Here we briefly present 


a three-state two-infection SI\SI 2 S model on arbitrary networks studied in [156] and further 


analyzed in 1571. 


The model is described as follows. Let Y/ be the number of neighbors of node i infected by 
the first disease I\. A node i in the susceptible state S transitions to the infected state I\ with 
rate f3}Y 2 . Similarly, a node i in the susceptible state transitions to the infected state J 2 with 
rate fifY 2 . Each node has its own recovery rate for each disease given by 5} and S'f. For example, 
a node i in the infectious state I \ recovers to the susceptible state at rate 5}. Figure [6] shows the 
interactions of this three-state SI 1 SI 2 S model. The deterministic version of this model is then 
given by 


N 


Pi = ~Pi ^iiPiPj 1 + PiPj*) + 5 iPi 1 + 5 iPi 2 ’ 


3 = 1 


N 


pI 1 = P S i a 'r j ip'/ - 
3 =1 

N 

2I 2 


pi 2 = Pi 


J2 ai jffi pI j 2 ~ s 'iPi 2 - 


3= 1 


For simplicity, we have only presented this model assuming both infections evolve over the 


same graph structure A. Instead, the works [156], 1571 provide analysis for these dynamics 
over possibly different structures. A few recent works have studied the problem of controlling 


multiple diseases in different scenarios 1571—[ 159 ]; however, these works are still in their infancy 
and there are still many open problems left to be solved. 


All the works about epidemics on networks we have discussed in this article have assumed a 
fixed graph structure. However, it is easy to surmise that this may not be a good assumption 
depending on the time-scale of a spreading process. For instance, in the context of diseases, 
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the network of contacts in a human population is constantly changing. Hence, a time-varying 
network model might be more appropriate, albeit more challenging to analyze. There is still 


very little work analyzing these types of time-varying models. Some recent works exist [160]- 


1 162) that have begun tackling these problems and laying the foundation for future works in this 
branch of epidemic research. As with optimal control, similar problems have been studied in 
different contexts such as information dissemination in mobile networks [163!], but far less has 
been considered in the context of epidemics thus far. 

In addition to the models we have presented throughout the article, it is worth mentioning that 
many works present the same types of models from a game-theoretic perspective [ |130[ , [ 1311, 


1164[-[166]. This is another space in which there is not yet a significant amount of study, but 


some seminal works have shown its usefulness in modeling spreading processes, especially in 


the context of control and optimization [ 167 ]—[ 170 ]. 


Conclusions 


This article has reviewed and analyzed some of the most popular models studied in 
epidemiology. In particular, we have presented deterministic and stochastic models in the context 
of both population and networked dynamics. We have described many results concerning the 
optimization and control of epidemic dynamics, while also outlining a number of new avenues for 
further exploration in this field. Although the focus of this article was on disease and epidemics, 
it should be emphasized that the same mathematical tools and results apply almost directly to 
a vast number of other spreading processes including information propagating through a social 
network, malware spreading in the World Wide Web, or viral marketing. 

Despite the vast literature studying the problems discussed in this article, there are many 
interesting control problems left to be solved, particularly those in the context of networked 
dynamics. There is plenty of work left to be done to really harness the power of these results 
and make a real societal impact; especially in understanding how to effectively control these 
processes on complex networks. In this respect, control engineers truly have a lot to offer in this 
reemerging field of research. 
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Infection rate 



Recovery rate 


Figure 1. Two-state Susceptible-Infected-Susceptible (SIS) model. An individual in the infected 
state / transitions to the healthy or susceptible state S with some recovery rate and from the 
susceptible state to the infected state with some infection rate. 
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Figure 2. Three-state Susceptible-Infected-Removed (SIR) model. An individual i in the 
susceptible state S can transition to the infected state / with some infection rate /3f ff and from 
the infected state / to the removed state R with some recovery rate . 
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Figure 3. Population dynamics of the two-state SIS model. These models assume a well-mixed 
population, meaning that each individual in the population is equally likely to contract a disease 
from anyone else in the population. An infected individual (red) naturally recovers at a rate 5 > 0, 
depicted by the red cross. A healthy individual (blue) if affected by each infected individual in 
the population with rate (3, depicted by the red arrows. 
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Figure 4. Network dynamics of the two-state SIS model. A node i has a natural recovery rate 
S, depicted by the red cross, at which it transitions from the infected state / to the susceptible 
state S and is affected by each infected neighbor j with rate /3, depicted by the red arrows. 
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5 


Figure 5. Three-state compartmental Susceptible-Protected-Infected-Susceptible (SPIS) model. 
An individual i in the infected state / transitions to the healthy or susceptible state S with a 
natural recovery rate 5. An individual in the susceptible state transitions to the protected state P 
at a rate fi(p s ,p ! ,p p ) that depends on the entire network state and to the infected state at a 
rate f3Yi proportional to the number of infected neighbors Y, that node i has. An individual in 
the protected state transitions to the infected state at a rate f3 0 Yi where /3 0 < /3 captures the 
fact that this individual is in a less susceptible state than normal, for instance due to behavioral 
changes or vaccination. 
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Figure 6. Three-state compartmental SI \SI 2 S model for two diseases. An individual i can be 
in the healthy or susceptible state S, or infected by one of two possible infections, but not both 
simultaneously. An individual in the first infectious state l\ recovers to the susceptible state at a 
natural recovery rate 5}. Similarly, an individual i in the second infectious state J 2 recovers at a 
natural recovery rate <5f. An individual in the susceptible state transitions to the infectious state/;,, 
at rate for k <E {1,2}, where Yf is the number of neighbors of i that are in infectious 

state //,. The parameter /3f captures the effect that neighbors of node i has on it for infection //,. 
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Sidebar 1 

Graph Theory 


A graph is a mathematical description of a given network. A graph consists of different 
nodes, or vertices, and links between the nodes, or edges, that describe the interactions between 
the nodes. In the context of epidemics, a the meaning of a single node depends on the granularity 
of the considered model. For example, a node at the lowest level can represent a single person 
and links to other nodes can represent the interactions this person has with others. On a much 
higher level, a single node can represent an entire city of people, and links to other nodes can 
represent the interactions this city has with others; for example, traffic flow between cities. See 
“Meta-Population Models” for further details. 

Formally, we define a directed graph Q = (V, E) as a pair consisting of a set of N 
vertices V and a set of edges E C V x V. The adjacency matrix A G M> 0 x7V of Q satisfies 
Oij = 1 if and only if (■ Vi,Vj ) G E. Edges are directed, meaning that they are traversable in one 
direction only. The sets of in-neighbors and out-neighbors of v G V are respectively 

A f in (v) = {v' G 1/ | (u» G E}, 

M™(v) = {v' G V | (v,v') G E}. 

We say that a graph is undirected if for all a l3 = 1, it is also true that a /t = 1. In this case the 
set of in-neighbors and out-neighbors for each node are identical. 

A directed path P, or in short path, is an ordered sequence of vertices such that any two 
consecutive vertices in P form an edge in E. A graph Q is strongly connected if for all vertices 
v G V, there exists a path to all other vertices v' G V. 
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Sidebar 2 

Meta-Population Models 


Throughout this article we often refer to individuals and discuss the state of all individuals 
in a network. However, especially in the context of diseases spreading through populations, the 
number of individuals N in a given network can be quite large. Instead of considering the 
entire population of interest together, meta-population allow groups of individuals to be lumped 
together into subpopulations under some assumptions. 

Consider the heterogeneous network SIS dynamics ([9]) 

N 

Pi = SiPi + ^2PjiPj(l - Pi)- (SI) 

3 = 1 

This was essentially introduced with pi referring to the probability that an individual i is infected 
(see “Networked Mean-Field Approximations” for further details). This means an A r -dimcnsional 
system must be analyzed to properly study how this model evolves which can be difficult if N 
is very large. 


Instead of studying the state of each individual in the population separately, we can 
instead create M « N subpopulations to approximate the dynamics of the entire N dimensional 
system. We can then keep track of the state of each subpopulation rather than the state of each 
individual in the population. This was originally done and analyzed for M — 2 and turned out 


to be easily extendable [171]. 


Let i e {1 denote the /th subpopulation with n, individuals, where each 

individual from the original population with N people is assigned to exactly one subpopulation. 
In other words, we have Y2jLi n j = N. Note that the number of individuals in each subpopulation 
do not need to be the same. 


We now define the dynamics of the meta-population model assuming that each subpop¬ 
ulation i is well-mixed and have homogeneous recovery rates <)[. In other words, within each 
subpopulation, each individual is assumed to have equal contact with everyone else. This is the 
same way the deterministic SIS population dynamics ([6]) are derived; however, we have the extra 
consideration that subpopulations can affect each other as well. The infection rate /i',- captures 
the effect that subpopulation j has on subpopulation i. Note that it is not required that /37 = f3[- 


51 



nor does it make sense to. Since subpopulations can have different numbers of people, it is 
reasonable to think that one subpopulation i can affect another subpopulation j more than j can 
affect i. Letting denote the fraction of individuals in subpopulation i that are infected, we can 
then write the dynamics as 

M 

Xi = -5[xi + y^ J Pj i x j ( 1 - Xi). (S2) 

3 = 1 


This has now reduced the original iV-dimensional system into an M-dimensional one. 
Furthermore, it might often make more sense to consider a meta-population model instead of an 


entire population one to begin with. To properly define the full network SIS dynamics ( |STj ), we 
require parameters that describe the natural recovery rates and interconnections of all individuals 
within the population. Instead, it is much more reasonable to believe these parameters can be 
estimated for entire groups of people and a reasonable meta-population model can be described 
with the same level of granularity. Furthermore, state information in the meta-population model 
can be estimated by looking at numbers of infected individuals in a given subpopulation compared 
to the total numbers of individuals rt,; in this subpopulation. For example, a node i € {1,..., M } 
at the lowest level of granularity simply recovers the full network SIS dynamics where each node 
represents a single person and links to other nodes represents the interactions this person has 
with others. On a much higher level, a single node can represent an entire city of people, and 
links to other nodes can represent the interactions this city has with others; for example, traffic 
flow between cities. 
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Sidebar 3 

Geometric Programming 


Let x e M> 0 , where .x'i,..., .x\y > 0 denote N decision variables. In the context of 
geometric programs, a monomial function h(x) is a real-valued function of the form //(x) = 
• • • x ]with c 0 > 0 and a* e M for alii e (1,..., iV}. A posynomial function g(x) is 
a real-valued function that is the sum of monomials, q(x) = c k x^ 1,k x 2 2,k . . . X a /’ k , where 

c fc > 0 and a iik e K. for alii e {1,..., N} and k e {1,..., Ji }. 

Before stating the definition of a geometric program, the following class of functions 
will be useful. 


Definition 17 A function / : M ;V —>• M. is convex in log-scale if the function 

F (x) = log/ (exp x), (SI) 

is convex in x (where expx indicates component-wise exponentiation). 


Remark 18 Note that posynomials (hence, also monomials) are convex in log-scale /172 ]. 


A geometric program (GP) is an optimization problem of the form 


minimize /(x) 

X 

such that qi(x) < 1, i — 1,..., m, (S2) 

hi (x) = 1, i = 


where / is a funcion that is convex in log-scale, g* are posynomial functions, and hi are monomial 


functions for i e (1,..., A r }. A comprehensive treatment of GPs is provided in [ 173 ]. A GP is 


a quasiconvex optimization problem [ ] 172 1 that can be transformed to a convex problem utilizing 
a logarithmic change of variables y, = log x t , and a logarithmic transformation of the objective 
and constraint functions. The GP in (|S2|) can then be written in the transformed coordinates by 


minimize F( y) 
y 

such that Qi (y) < 0, i — 1,..., m, (S3) 

bjy + log di = 0, i = 1,... ,p, 

where Qi (y) = log g*(exp y) and F (y) = log/ (exp y). Also, given that h % (x) = 

diX^x^' 1 ... x^'\ we obtain the equality constraint above, where fe* = (bij ,..., 6jv,*)- 
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Since / (x) is convex in log-scale, F (y) is a convex function. Furthermore, since q t 
is a posynomial (and therefore convex in log-scale), Q l is also a convex function. This shows 
that (fS3T) is a convex optimization problem in standard form and can be efficiently solved in 


polynomial time [ 1721. 
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Sidebar 4 

Networked Mean-Field Approximations 


The method of going from a stochastic model to a deterministic mean-field approximation 
is certainly not one that can be overlooked. The derivations of these approximations, their 
accuracy, and what they can tell us about the original stochastic models is an entire area of 
research all by itself. 

We discuss briefly here how one can go from the stochastic model ([7]) to the deterministic 
one ([8]). Recall the stochastic model 

Xi : 0 ->• 1 with rate fd X v 

X; : 1 —>• 0 with rate 5. 

Given the entire state X(t) at some time t, we can write the probability of state i at a future 
time t' — t + At for small At by 


P{Xi{i!) = 0| Xi(t) = 1 ,X(t)) = 5 At + o(At), 

P{Xi{lf) = l\Xi(t) = 1 ,X{t)) = l-SAt + o(At), 

P(Xi(lf) = l\Xi{t) = 0 ,X(t)) =pY. X ^) Xt + °( At )> 

jeXi 

P(Xi(t') = 0|7Q(f) = 0,X(t)) = 1 X^t)At + o(At). 

j&Xi 

As we take At to 0 in these forward Kolmogorov equations, we can write the exact dynamics 
of the expectation by 


dE\Xj) 

dt 


-E 


s + ( l-xjp^xM 

jeATi 


-E[Xi]8 + E 


(1 - X,)P v x s (t)) 

jeXi 


The complication now comes from the term E[X t Xj\ relating the covariance of the random 
variables X t and X :) with their independent probabilities. The mean-field approximation ([8]) (and 
similar ones for different variations of the stochastic model) are then obtained by assuming that 
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E[XiXj\ = E[Xj]E[Xj\ for all i ^ j. In other words, we are assuming all the random variables 
are independent. 


However, this is clearly not true. This means that for any fixed population with a stochastic 
model, the deterministic approximations we study are just that; approximations. Naturally, this 
begs the questions of how accurate they are in describing their stochastic counterparts. 


Although we only consider an approximation of the expected values it has actually 


been shown that these are upper-bounds on the actual probabilities [72], [73], [781 (this is 
essentially done by showing that E[XiXj\ > 0 for all i ^ j). Fortunately, this has very positive 
implications on attempting to control the underlying stochastic process by using the deterministic 
mean-field model: By stabilizing the deterministic approximations, we can make claims like the 
ones presented in Remark |7J More specifically, if we can ensure that the disease-free equilibrium 
of the deterministic model is globally asymptotically stable, then the stochastic system will reach 
the disease-free absorbing state in sublinear time (with respect to the size of the network) in 
expectation. 


In [174], the authors begin looking at how accurate the deterministic mean-field 
approximations are in describing the stochastic models, rather than just guaranteeing the upper 
bound. This is still an open problem for arbitrary networks. 


work 

more 


Furthermore, all works above only consider the SIS dynamics. 


[175] provides this type of analysis for a three-state SIRS model, 
complicated models in general are still unsolved problems. 


Although the recent 
rigorous analysis for 
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